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Abstract. This paper introduces a new class of Cox models for dependent bivariate data. 
The impact of the covariate on the dependence of the variables is captured through the 
modification of their copula. Various classes of well known copulas are stable under the 
model (archimedean type and extreme value copulas), meaning that the role of the covariate 
acts in a simple and explicit way on the copula in the class; specific parametric classes are 
considered. 

I—* 

1. Introduction 

m 

The aim of this paper is to present a new description for bivariate dependence. It extends 
f— i \ the proportional hazard (PH) model and is relevant in various fields of biostatisctics and 
industry. Denote z an environmental covariate and consider two positive random variables 
X and Y with absolutely continuous survival functions 

(1.1) F\x) = ¥(X > x; z) 

and 



c3 



i/V (1.2) G (y) = P(y > y; z) 

under the covariate z. The joint distribution of (X, Y) is modelled through (11 .ip and (jl.2p 
and through the conditional s.d.f 



^1- 



(1.3) Gl(y) = ¥(Y>y\X>x;z). 

The precise setting of our model is as follows: we assume that the survival function of X 
depends on z through a PH model, i.e 

(1.4) X x (x) := ~ \ogF\x) = \° x (x)$(z) 

\AjJu 

for some positive function z —*t &(z) and some baseline hazard \*x(x) corresponding to z — 0. 
Denote 

Xy\x>x(v) ■= ~^ log ^^)- 
The conditional survival function of Y depends on z through the PH model 

(1.5) ^y\x>M = >?Y\x>Ay)y 

for some positive function z — > ^{z) and some baseline hazard Ayi x>x . Denote 



1.6) (Ml) 



X z x (x) = \° x (x)${z) 
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Other models for b i yariate depe n dence hav e been d e fined, for example (jClaytonl . 11978 



Clayton and Cuzickl . Il985l ; lOakesl . Il989l ). In IClaytonl (ll978[ ) where joint lifetimes of sons 



A 


[y 


X = x) 


A 


(y 


X > x) 



and fathers are considered, a Cox model is used to handle the role of the covariate upon the 
margins, and the association between the margins is described through the fact that the local 
association measure 

6 (x,y) :-- 

is independent upon x and y and upon z. This index expresses the influence of parental history 
of a given disease upon the incidence in the offspring. X is the lifetime of the father, while 
Y describes that of the son. It is assumed that association arises because the two members 
of a given pair share some common influence and not because one lifetime influences the 
other. In fact these models are of frailty type, and the margins are independent conditionally 
upon the covariate. More globally the standard model building strategy is to have m a rgina l 



su rvival funct i ons a nd a copula for their depend ency; see f.i. IChoi and Matthews! (l2005[ ) 



or 



Saidetal 



( 120091 ). In Lindstrom et al.l ( 2007 ) a study on the familial concordance in 
cancer survival based on a Swedish population showed that cancer specific survival in parents 
predicts survival for the same cancer in their children. The risk in dying in children in 
relation to parental survival was modelled by use of two PH models; first parental survival 
was modelled and next survival risk in children in relation to parental survival was assessed; 
thus in this case the lifetime of the father influences the lifetime of the son whatever the level 
of the covariate (here z represents the different cancer sites), however in relation with it. In 
this case Model ( II. 6p seems to be natural. Setting x = in (II. 5p shows that Y follows a 
PH model so that (II. 6p is PH on both components. The copula of the couple of r.v's (X, Y) 
for a given z can be written in terms of the values of $(-z) and ^{z), as seen further. The 
unusual feature of this model is positive in a number of cases, since z acts on the margins, 
but also specifically on their dependence. The claim that the goal of copula modelling is to 
distinguish the parameters for the dependency from those associated to the marginal models 
cannot be considered as a general principle. Other authors have considered cases when th e 



association of the margins is specifically related to their distributions; see f.i. iGuptal (120081 ) . 
In conjunction with (jl.6p it is worth noting that direct approaches based on regression type 
models cannot satisfy our purpose. Indeed consider for example a model defined through 

" X = r(z,U) 
Y = s{z,V) 

with r(z, .) and s(z, .) strictly increasing for all z. Then following iNelsenl ( 120061 ). Theorem 
2.4.3, (X, Y) has the same copula as (U, V) for all z, which implies that the covariate z plays 
no role in the dependency of X and Y. Introducing a new model imposes to determine its 
range of applicability; it will be shown that (jl.6p is adapted for positive quadrant depen- 
dence (PQD) between the margins, a concept which is recalled in Section [2j Stating that 
the margins follow a PH (typically Cox) model can be checked using standard tools (see 
f.i. Grambsch and Therneauj. ll994T) . PQD property can be tested through a Kolmogorov- 



Smirnov test (see f.i. IScailletl . 120051 ). Let us now show the main results which we present in 



connection with model ( II. 6ft : 

(1) The TP 2 class of sdf's is a subclass of the Positive Quadrant Dependence (PQD) sdf 
class and is stable under the model which is properly defined when the hazard baseline 
H is TP 2 . This class appears quite naturally as the one under which the model is 
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(2) 



properly defined, and it is appropriate for the description of positive dependence 
between its margins. Definitions of TP 2 and PQD properties are given in Section [2J 
Since the TP 2 property of a multivariate sdf refers only to its copula, model (II. 6p 
describes the changes of the baseline copula induced by the covariate. Also this 
implies that the model is valid independently of marginal distributions. Only the 
structure of dependence is involved in the domain of validity of the model. 
(3) Two main classes of copulas are stable under the model namely: when the baseline 
bivariate copula is in such a class, so is the copula for all value of the covariate z. 
The class of extreme values copulas (eve) enjoys this property. The class of extended 
archimedean copulas is also stable under the model. This class results as a special 
by-pro duct of a techniqu e inten ded t o produce asymm etric copulas due to Genest et 



al.; see iFrees and Valdezl (120071) and iLiebscherl (l2008h . The so-called class of logistic 



asymmetric copulas (see iTawnl . Il988l ) . which is a simple extension of the Gumbel 
family of copulas, enjoys an important role in the present model. It is stable under 
the model and admits a simple parametrization. The covariate z acts in an adaptive 
way when the value of the covariate is changed. It is the only bivariate distribution 
in the class of frailty models which enjoys such properties in the model. 

This paper is organized as follows. In Section [2] we briefly recall the necessary background 
from bivariate dependence. Section [3] describes the model. In Section 0] we focus on the 
asymmetric Gumbel class of copulas, which is the natural parametric setting of our model; 
we also provide some connection with bivariate frailty models. All proofs are deferred to the 
Appendix. 



2. Some useful facts in bivariate dependence 

Let X and Y be two random variables (r.v) with joint sdf H, with margins F and G. All 
dependence properties of X and Y are captured through the survival copula C which is a cdf 
defined on [0, 1] x [0, 1] through 

C(u,v) = H(F^{u),G^(v)) 

where u and v belong to [0, 1] and where F (t) : = sup{x : F(x) > t\. I t is eas ily checked 
that C is indeed a copula. The definition of a copula is given in iNelsenl (120061 ). definition 
2.2.2. We will make use of the following definition and notation. 

Definition 1 (Min-id property). A bivariate cdf H is min-infinitely divisible (min-id) if for 
all positive 7, IP is a sdf. 

Assume that H is min-id and let V = (X, Y) be a random vector with sdf H. Then for all n 
in N, H is a sdf. Further let (Xj, Yi) , i — 1, . . . , n be n copies independent and identically 

1 In 

distributed with sdf H . It holds 

V= (min Aj,min Yj). 

i i 

Definition 2 (PQD property). X and Y are positively quadrant dependent (PQD) iff, for 
all (x,y) in IR 2 ,P(X > x, Y > y) > F(X > x)F(Y > y); in this case we also say that H is 
PQD. 
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Definition 3 (TP 2 property). A mapping from R onto R is totally positive of order 2 



™2 



and 



(TP 2 ) if (f)(x,y) > for all (x,y) in 
4>{xi, V2)4>(x2, Vi) > 0, for all x x < x 2 and y x < y 2 . 
Remark 1. When <p is C 2 , then <p is TP 2 iff 



'x 2 ,yi) 0(^2,2/2) 



<f>(xi,yi)(f>{x2,y2) - 



(2.7) 



dx ^ dy ^ dxdy 



x,y)<l){x,y) 



The proof is given in Resnickl ( 1987 ). p. 254. 
We also recall the following results. 



Theorem 1 ((Jod (119971 ). Theorem 2.3). If H is a TP 2 sdf then H is PQD. 



Theorem 2 (jJod ( 119971 ). Theorem 2.6). Let H be a cdf, then H is min-id iff H is TP 2 



The relation between properties of the s.d.f 's and their copulas is captured in the following 
result. 

Lemma 1. Let H be a sdf with copula C. Then H is TP 2 iff C is TP 2 . 

Proof. By Remark [H H TP 2 =>- ^(x,y)^-(x,y) < §^(x,y) x H(x,y). Furthermore, 
H(x,y) = C(F(x),G(y)). Therefore, 



dH ( 

dx 



x,y) 



-—(u,G(y)) 

OU u=F{x) 



dH 

-dy~^ v) 



v=G(y) 



x g(y) 



d 2 H 

dxdy 



x,y) 



d 2 C 
dudv 



u= f [x) x f(x) x g(y) 

v=G(y) 



Hence, 



du 



u,G(y)) 



u=F(x) dv 



x -^(F(x),v) 



v=G{y) dudv 



u=3x) xC(F(x),G(y)) 

v=G{y) 



Definition 4 (Archimedean copula). An Archimedean copula is a function C from [0, l] 2 to 
[0, 1] given by C(u, v) = <y?' _1 l ((p(u) + (p(v)) , where ip is a continuous strictly decreasing convex 
function from [0, 1] to [0, 00] such that <p(l) = 0, and where yj' -1 ' denotes the "pseudo-inverse " 
of if, namely 

J-i]u)-S ^\t) fortin[0,<p(0)] 
* W ~\ fort><p(0) 

When (f(0) = 00, (p is said to be strict and (p}~^ = tp^ 1 . <p is called a generator. 



The class of so called extreme value copulas (eve) is important in this model, although not 
related here with the theory of bivariate extremes; therefore we define an extreme value 
copula through the basic Pickands representation, without further reference to the theory of 
bivariate extremes. 
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Theorem 3 (Pickands Theorem). C is an extreme value copula iff there exists a convex 
function A defined on [0, 1], which satisfies A(0) = A(l) = 1 and max(t, 1 — t) < A(t) < 1 
such that 

logv 



(2.8) C(u,v) = exp log(uv)A[ 



Jog uv / 

The function A is referred to as the dependence function or Pickands function of the copula 
C. 

3. Introducing covariates in dependence models 

3.1. Description of the model. Not all baseline survival d.f's H° defines a model, so that 
X z x and Xy\x>x are the marginal and conditional specific cause hazards for some bivariate sdf 
H z with margins F* and G* under a given covariate z. We conclude from the first equation 

of ( II. (jp that F (x) = [F (x)\ . By the second equation in (11.6j) . plugging x = 0, we 

get G z (y) = (^G°(y) Ss j . The model is defined when z holds if (^F°(x)J {H Ylx>x (y))^^ 
defines a sdf. Notice that 



(3.9) = (#W)) " (F°(x) 

which is indeed a sdf when $(-2) > ^f(z) > and ( H [x, y)J is a sdf. 

Also not all bivariate survival d.f's H° are such that for all positive 7, (^H°^j is a sdf- Min- 
infinite divisibility of the baseline hazard seems to be a natural assumption here. Assume 
therefore that: 

(H) H° is min-infinitely divisible 

By Theorem [2] and Lemma (TJ (H) holds iff C-^o is TP2. We have the following result. 

Proposition 1. When fH\) holds then H* defined in ( Iff. ,91) is a sdf for all z such that <3>(z) > 
> 0. 

Let us consider the case when < §(z) < ^(z). Analogously with ( II. 6p the model may then 
be written 



(3.10) (M2) 



/ K(v) = KivMz) 

S \z f „\ \o f„\f 



permuting the role of X and Y. In a similar way to the above we have that 

#>,y) = (V (*,</)) (G°(y)) 
is a proper sdf. To summarize the above arguments we state: 

Let the model be defined by flUE) if > *(z) and by (l3~T0|) if $(z) < Call 

(M) the model defined through 

(M) := (Afl)l« W >* W + (M2)l* w< * w . 



This model is well defined, even if *&(z) and ^(z) are not ordered uniformly on the covariate 
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z (which can be multivariate); the functions <l> and \1/ can be easily estimated through the 
data, since they characterize the marginal PH models in (M). Suppose that X and Y are 
fitted to the same scale under the baseline, namely F (t) = G (t) for all t. Then > ^f(z) 
implies F (t) < G (t) for all t, stretching the fact that X becomes stochastically smaller 
than Y under the stress parameter z. 

Identifiability of (M) holds; assume for example that Q(z) > ^f(z). Then &(z) and ^f(z) 
are defined uniquely. Indeed, assume 

H z (x,y) = (jF(s,y)) 



(lf(x,y)) (Z> (f\x 



for all x, y. Then taking logarithms yields $(2) = $'(z) and ^(z) = ty'(z). 
When flHl) holds then for all z , H is a sdf and 



(3.11) H z (x,y) = l Hz) > Hz) (lf(x,y)y iZ) ^F\x)y {Z) ^ 

+!*(*)<*(*) l-ff (x,y) J [G (y) 



Min-infinite divisibility of the baseline will also make any H z min-infinitely divisible, showing 
that this class is stable under (M). Indeed for any positive 7, 

(tf(a;,?/)) 7 = l* w >* w (tf (x,y)J (F (x)J 

/ o \7*W /_o \7[*(*)-*(*)] 

+ 1 0(z)<\E r (z) 

which still is a sdf. By Theorem [2] and Lemma [TJ min-infinite divisibility is not a property of 
the cdf but of its copula. Formula (13. lip can be written for copulas through 

$(z)-*(z) / _J_ 1 \*M 

(3.12) Cjj*(u,v) = l* w >*(*)« *w C^o ?;*<*> J 

(jH|) is only a sufficient condition for the model to be defined. The following example illustrates 
this fact. 

Example 1. Let C-go(u,v) = w exp{- fllogitlogu}, with 6 £ (0;1]. This is the Gumbel- 



Barnett family (see Nelsen , 200d . p. 119). By \2. 7]) it is easy to check that CLjo is not TP2. 



Using $3.12\) and assuming > ^f(z) we obtain 

Cjjz(u,v) = uvexp | — ^ ^ log u log v | 
which still is a Gumbel-Barnett copula when belongs to (0; 1]. 

This example shows that (|H|) is indeed the only acceptable condition for existence. Otherwise 
the baseline hazard H defines a model only for specific values of the covariate. This motivates 
our interest in good classes of min-infinitely divisible copulas which we intend to regress on 
the covariate z. 
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3.2. Stability properties of the model. We introduce two classes of copulas which are 
stable under (M). 

3.2.1. Extended archimedean copulas. Among all possible types of bivariate dependence which 
can be described through the present bivariate Cox model, there exists a class of copulas which 
contains the archimedean copulas and which enjoys peculiar stability properties. 

From now on, denote 

(3.13) a(z)=min -f4,l 

\§{z) 

and 

(3.14) P(z) = min ( 1 



Let ip be a generator and C vo be the archimedean copula with generator (p . We assume 
that C Vo is TP 2 . When H is defined, set C-jp its copula. 

Proposition 2. Let H° be a sdf with TP 2 copula C vo . Then, H* is defined for all z in the 
domain of $ and ^ . Further 

(3.15) Ctj*(u,v) = u^v^cp; 1 (<p M (u^)+<p M (v^)) 
with 

(3.16) <p g (t) = v? ^*w«wj = Lp {t^m^j . 

More generally we have, denoting IL(u,v) = uv 
Proposition 3. Assume that 

(3.17) C^iu.v) = Il{u l - K ,v l - r ')C^(u K ,v n ), < «, rj < 1 
£/ien 

(3. 18) Cg- (u, u) = Uiu 1 - ^, v 1 -^)^ (u a{z)K , v 0{z)v ) 

where C Vz denotes the archimedean copula with generator ip z defined in $3.16}) . 

3.2.2. Extreme values copulas. We show that the class of eve's als o enioys stability p roperties, 



6.Z.Z. Extreme values copulas, we snow mat tne class oi eve s als o enioys stability p i 
as seen in the present Section. Extreme values copulas are TP 2 ( Hurlimannl . 20031 ) . 



Proposition 4. When H has an eve CW) with Pickands function A then, denoting Cjj" the 
copula of H z and using f)3.12p . we have 

logf 



(3.19) C-jjz(u,v) = exp 



\og(uv)B z 



log uv 



with 
(3.20) 

B z (s) = 1 - W(z)K(z) - sW(z)[l - K(z)} + W{z)[{\ - s)K(z) + s}A 
where K(z) = andW(z) = min f^y, l) = P(z). 



K(z)(l-s) + s 
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Remark 2. This basic result shows that C-^ is an eve with Pickands function B z . Proposi- 
tion^ shows that the class of eve's is stable under (M). Although the copula of H is an eve, 
this does not imply in any respect that its marginals should be extreme value sdf's. 

From f)3.20p we deduce the transition formula which links B z to B z for two different values 
of the covariate. It holds 

Proposition 5. Under (M) let H° has an eve. Then with the above notation, for all z, z' , 
(3.21) 



a(z) \ct(z) (3(z) 



a(z') J(z') 
v J a(z) + I3{z) 



6 m 



' a(z) p(z) 



Proposition proves that the transition from z to z' is independent of the baseline depen- 
dence function. Formula (I3.2ip can be seen as a kind of expression of the proportional hazard 
property, which links two hazard rates independently on the baseline. When the covariate 
acts equally on X and Y, i.e. &{z) = ^{z) for all z, then B z (s) = A(s) for all values of s as 
seen in Proposition HI Thus, the copula of H equals that of the baseline H ; the dependency 
structure of X and Y should not be altered through (M). Only the marginal distributions of 
X and Y in this case reflect the role of the covariate. 



We propose some illustration. We use $(2) = e az , with a = 1.5 and ^(z) = e' 32 , with 
(3 = 2. Figure [37T1 illustrates formula f )3.15p . We represent the change of the density of 
Cjjz with z. The archimedean copula is the Clayton copula whose generator is defined by 
(p(t) = t~ 9 — 1. We take 6 = 3. In this figure the model tends rapidely to independence 
since the density of the copula tends to 1 as z increases. Figure 13.21 illustrates the transition 
formula f)3.20p . The baseline copula is the Gumbel copula with 6 = 3. The Pickands function 

of the Gumbel copula is A(t) = \t e + (1 — . The dependence functions are ordered wrt 
z. As z increases, the model tends to independent marginals. 



9 





FIGURE 3.1. Illustration of (13. 18[) with the Clayton copula density for the 
baseline (z = 0) 



4. Asymmetric logistic models of dependence 

This section deals with specific parametric models for dependence which are stable under (M). 
We consider model (M) specialized in the case when the copula of if is a Gumbel copula. 
The margins of H can be any. It is a simple parametrized model of copulas, which is an eve 
on one hand, and which models frailty bivariate dependence, being hence a n archimedea n 
copula. Indeed it is the only copula satisf ying jointly these two properties (see iNelsenl (120061 ) . 
Theorem 4.5.2; iGenest and Rivestl ( 119891 ). statement A). 
The Gumbel copula writes 

j(-logw) e + (-logwfj 1 



C(u, -y):=exp 



with 9 > 1. The dependence function of this copula is 



(4.22) 



A(s)= s e + (l-sf 



1/0 



Assume that H has an eve with dependence function A. When the covariate z acts, the 
dependence function B z defined through Proposition [5] determines the asymmetric logistic 
copula. This copula has three parameters a(z),(3(z) and 9. Recall from ( 13. 13ft . (13. 14f) and 
proposition 0] that 



no !U)!1 l^y,l 
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Figure 3.2. Illustration of formula (I3.20|) 



and 

It holds 
(4.23) 



mm 



B z {s) = 1 - a{z) + [a(z) - 0{z)\ s + a{z) e (1 - s) 6 + 0{z) 



When a(z) = (3(z) (which implies that they equal 1), i.e. when z acts equally on X and Y, 
B z (s) = A(s) for all s. The copula of H is 

(4.24) Ch*( u ,v) = ni^-^y-^ciu^y^). 

where H(u, v) := uv is the product copula (see Nelsenl . 2006 . p. 11). The dependence function 
B z is an asymmetric form of t he Gumbel d ependence function A defined in ()4.22|) . This is the 
asymmetric logistic model in Tawn ( 19881) wh en the margins are standar d exponential. As 
develo pped by Khoudraji (see lKhoudrajil . ll995l . chap 4) and Genest et al. in lFrees and Valdez 



fl2007f ). Proposition 3, given two dependence functions A\ and A 2 , two constants k and 7] 
with < K, 7] < 1 the function defined through 



B(s) := (ks + r)s) A 1 



KS 



+ (ks + rjs) A 2 



KS 



KS + T]S J \KS + TjS / 

where s denotes 1 — s, is the dependence function of the extreme value copula defined by 



C 



aAu 



l-K 



)CA 2 (u K ,v' n ). Genest et al. define this procedure as a technique to generate 
asymmetric copulas. The class of copulas defined in (|4.24p has been introduced by Genest 
et al (see their Proposition 2 in Frees and Valdez ( 20071 )). 
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We now analyze this class of copulas in terms of frailty models. The Gumbel copula is 
associated with a frailty model of order 1, namely 

C(u,v) = A -1 (A (it) + A(u)) 

where A" 1 ^) = e- sw dM 1/e (w) = e~ s 1 , 9 > 1, is the Laplace transform of the positive 
stable law My$ on IR + with tail heaviness index 1/9, lo c ation parameter 0, scale par ameter 



1 and skewness parameter (see Ravishanker and Dey ( 2000[ ) and the example 5 in Oakes 
( 1989[ )). Denote W a positive random variable with cdf Myg. A bivariate sdf H° with Gumbel 
copula C writes 



H°(x,y) 



{F(x)G(y)} w dM_ 



l/e 



w 



for some sdf F and G. Therefore if is a frailty bivariate sdf, with stable frailty measure and 
margins f °° F(x) dMx/e(w) and J °° G(y) dM\/e{w). 

Let U\ and U 2 be two independent r.v's, both independent of W. We assume that U\ and U2 
have a positive stable law M 1 and M 2 with tail heaviness index 1/9. The r.v U\ (resp. U2) 
has shape parameter — 1 (resp. — 1). Both have location and skewness parameters 

0. Define Si := Ui + W, i = 1, 2 (see lDrouet and Monbetl . |2004 - Denote ^{s) the Laplace 



transforms of the distribution of Si. Denote further ?/'~ 1 (s) the Laplace transform of the 
distribution of Ui. Let M denote the probability measure of (Si , 82). For arbitrary sdf H\ 
and H2 define the bivariate sdf 

H(x,y):= J J {H^x)} 31 {H 2 (y)} S2 dM(a u a 2 ) 

which we call a frailty model of order 2 si nce it implies a bivariate lat ent variable. Frailty mod- 
els of order two have been considered in Marshall and Olkinl (Il988h (see their formula (2.2)). 
The marginals of H are F(x) = J °° (ifi(x)} 1 dMs^si) and G(y) = / °° {H 2 (y)} 2 <iMs 2 (s 2 )- 
We prove that the copula of H is (I4.24p . Indeed 



H(x,y) = J {H l (x)} Ul dM\ Ul ) j {H 2 (y)} U2 dM 2 (u 2 ) J {H l (x)H 2 (y)} w dM x 



[w 



Introducing the Laplace transforms defined above and rewriting the marginals F(x) = 
ifi 1 (— logifi(x)) and G(y) = (f^ 1 (— log H 2 (y)) we obtain the following expression for the 
copula of H 



(4.25) 



C(U,V) = (^(u)) ^ {(pz{v)) A' 1 (<P!(U) + ip2{v)) 



since H(x,y) = C(F(x),G(y)). Substituting ^ and (p i7 i = 1,2 by their expressions in the 
above expression, noting that ip^ 1 = ■0r 1 A -1 , (I4.25P coincides with (l4.24j) . We now prove 
that for an adequate choice of Hi and H 2 the bivariate sdf H has same marginals as H . 
Indeed let 



H x {x) 
H 2 (y) 



cxp 



exp 



- - min ($(z), ^f(z)) logF (x) 



- ( -min (®(z),^(z)) log G°(y) 
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which yields F(t) = F (t) and G(t) = G (t) for all t. Therefore H and H coincide. We have 
proved 

Proposition 6. When H° is a frailty bivariate sdf with Gumbel copula, then for all z, H z is 
a frailty sdf of order 2 with asymmetric logistic copula given in {1^.2$ . 

More generally we have 

Proposition 7. The class of second order frailty models with asymmetric logistic copula is 
stable under (M). 

Proof. Let A denote the dependence function of an asymmetric logistic copula. 



A(s) := 1 - « + (« - 77) s + k° (1 - s) e + n e s 



l/e 



By (jHZDD it holds 



with 



B z {s) = 1 - k' + («' - rf) s + (1 - S y + rf 



Ot(z)K 



J9 J 



i/e 



These new parameters are in (0, 1) , cts Kj and rj. We have proved that the class of sdf with 
asymmetric logistic copula is stable under (M). Any sdf with such a copula is necessarily a 
frailty sdf of order 2. Indeed this follows from (I4.25P which enables identifying the frailty 
measure M of H as the joint distribution of (5i, S2) as defined here above. ■ 



Remark 3. It can be seen that the only sdf which are frailty of order 2 with eve are precisely 
the frailty models with asymmetric logistic copula C(u,v) = U(u 1 ~ K ,v 1 ~ v )Cg(u K ,v v ) with 
< r], K < 1, and where Cg is the standard Gumbel copula with parameter 9 > 1. 

5. Simulation results 

5.1. Stability of the estimate under the model. For a given copula Cg in a parametric 
family we simulated iV independent couples i = 1,...,N, with joint distribution 

function Cg. Estimation of 9 was performed using plug-in technique, leading 9. We repeated 
the procedure 1000 times. We compared Cg with C| for various z, given known functions $ 

and Here $(2) = e az with a = 1.5 and ^(z) = e l3z with (3 = 2. The figures hereunder 



show the mean relative error J J 



C$(u,v)-CI(u,v) 



dCg(u,v) with respect to z, together with 



Cl(u,v) 

the 95% confidence interval. This indicator is obtained on a grid of 100 points in [0, 1] x [0, 1]. 
Figure [5731 pertains to the Clayton copula Cg(u,v) = (u~ e + v~ e — 1)~« with 9 = 3. The 
estimate 9 is defined through 9 = where 

. . „ number of concordant pairs — number of discordant pairs 

5.26 t = t 

/ N 

2 



is the empirical estimate of the Kendall's tau (see iNelsenl . 12006k p. 158). Figure 15741 pertains to 
the Gumbel copula Cg(u,v) = exp[— {(— logu) e + (— log v with 9 = 3. The estimate 9 is 
defined through 9 = j^. It appears from those curves that a good estimate in the reference 



BIVARIATE COX MODEL AND COPULAS 



13 



zone (z = 0) propagates accordingly to other zones (indexed by z), without deteriorating the 
estimation accuracy. These facts also hold for very small values of iV; obviously the larger 
N, the better the accuracy, for all z. 

Remark 4. 

(1) These simulations are closely related to some industrial application where the envi- 
ronmental variable z has range [0,0.3]. Obviously a (resp. (3) is fitted accordingly 
around 1.5 (resp. 2), as in these simulations. Changing the scale of z modified the 
value of a and (3. This entails that no comparison can be performed with higher values 
of z on these graphs, keeping in mind the underlying applied statistical context. 

(2) In these simulations the paramater 6 of the copula is estimated using the baseline 
group. The estimation of the copula under z is obtained through a transformation 
of the copula under z — 0. Taking into account all the pooled data for all z values 
introduces complex estimation procedures. 
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Figure 5.3. Relative error pertaining to Clayton copula 
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N = 100 
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N = 1000 



N = 10000 




o 




Figure 5.4. Relative error pertaining to Gumbel copula 
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5.2. A case study through simulation. In this section, z is bivariate, z = (zi,z 2 ) with 
Zi G {0; 1}, i = 1,2. We consider model (M) with Weibull marginals. In the reference 
zone, z' = (0,0), the marginals are W{2, 12000) for X and 1^(1.5,8000) for Y. The copula 
in the reference zone is Clayton (archimedean) or Gumbel (which is both archimedean and 
e.v.c). The function $ (resp \1/) is exp(o: / z) (resp exp(/3'z)) where ex! = (0.1,0.06) (resp f3' = 
(0.07,0.25)). We simulated N = 200 couples with Gumbel copula or with Clayton copula 
under z' = (0,0), 100 values under z' = (1,0) an d 100 values unde r z' = (0, 1). Simulations 



was performed using Khoudraji's algorithm (see iKhoudrajil . Il995h . We estimated a and /3 
using the whole sample (400 couples) through Cox's partial likelihood estimation method 
(implemented in R via the coxph procedure); therefore we estimated ex and (3 independently 
by working separately with the XJs and the Y^'s. The copula parameter 9 under z' = (0, 0) 
was estimated through the plug-in of the empirical Kendall's tau. We used formula ( I3.15P or 
formula (I4.24p in order to obtain an estimated copula C z under z' ^ (0, 0) (with 3 parameters 
instead of 1) by replacing a(z), /3(z) in (13 . 13|) and (I3.14p and 9 by their estimators; the 
estimator of 9 is unchanged and has been obtained under z' = (0,0). On the other hand, 
these same formulas provide the theoretical copula C z under z. We performed the relative 

accuracy of the estimation scheme through J J — ^fej^r^^ dC z (u,v). The procedure is 
repeated 1000 times. The results are given in Table [TJ We give the average of the uniform 
relative error over the 1000 simulations (bold character) and its 95-percent confidence interval. 





z' = (0,0) 


z' = (l,0) 


z' = (0, 1) 


Clayton (9 = 3) 


0.93% 

[0.88%,0.97%] 


2.03% 

[1.95%,2.12%] 


2.60% 

[2.49%,2.71%] 


Gumbel (9 = 3) 


0.81% 

[0.77%,0.85%] 


1.55% 

[1.4996,1.62%] 


2.02% 

[1.93%,2.11%] 



TABLE 1. Relative error for different copulas 





z' = (0,0) 


z' = (l,0) 


z' = (0,l) 


Clayton 


3.49% 

[3.32%,3.66%] 


5.73% 

[5.44%,6.02%] 


11.65% 

[11.08%,12.23%] 


Gumbel 


2.35% 

[2.24%,2.46%] 


4.01% 

[3.82%,4.22%] 


8.50% 

[8.11%,8.89%] 



TABLE 2. Relative error on Spearman's rho 
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The estimators of $ and \I/ however are of mean accuracy, which nevertheless does not 
deteriorate the quality of the estimators of the copula, in all cases which we considered. 
In Table Q] the indicator mixes both propagation and estimation errors. The relative error 
w.r.t the theoretical propagated copula is very small in the case of the most commonly used 
baseline copulas. Measured in terms of Spearman's rho, the relative error does not deteriorate 
significantly either, as seen in Table [2j 

5.3. Propagation of misspecification errors. 

5.3.1. Assuming $ and \1/ known. We simulated N = 200 couples of r.v's with Clayton 
distribution function Cq with parameter 6 = 3. Here $ and \l/ are as in Section 15.11 We 
estimated the Kendall's tau through the classical non parametric estimate (I5.26p . We used 
a misspecifi ed mod e l assu ming that the data have been generated under a Ali-Mikh ail-Haq 
copula Gq { Nelsenl . 2006 . Table 4.1 p. 116) with parameter 9 = As seen in Nelsenl 

fl2006h . both copulas present common features. Misspecification may therefore occur as 
presented here. For various z we used formula ( I3.15P to define both the true copula C§ 
and the misspecified estimated copula G%. The misspecification error is defined through 



Err(z) 



II 



C*Ju,v) 




I 1 1 1 1 1 1 

0.00 0.05 0.10 0.15 0.20 0.25 0.30 



z 



FIGURE 5.5. Relative error pertaining to the misspecification when $ and \I/ 
are known 
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5.3.2. Assuming $ and ^ unknown. We simulated N = 200 couples of r.v's with Clayton 
distribution function C$ with parameter 6 = 3. z is bivariate and the functions $ and \I/ 
are as in Section 15.21 We estimated $ and \I> through partial likelihood and we used a 
misspecified model assuming that the data have been generated under a Ali-Mikhail-Haq 
copula Gg. The results in Table [3] show that the misspecification error is of order 17% and 
keeps stable through the propagation. At the contrary Table [1] shows that the error under 
the true model is much smaller and propagates with great accuracy. This enlights the need 
for a good specification in this model. 





z' = (0,0) 


z' = (l,0) 


z' = (0, 1) 


Clayton - AMH 


17.38% 

[17.36%,17.40%] 


17.61% 

[17.54%,17.67%] 


14.72% 

[14.65%,14.80%] 



TABLE 3. Relative error Err(z) pertaining to the misspecification when $ and 
\l/ are unknown 





z' = (0,0) 


* = (1,0) 


z' = (0, 1) 


Clayton - AMH 


53.39% 

[53.31%,53.46%] 


54.85% 

[54.61%,55.09%] 


55.29% 

[54.83%,55.75%] 



TABLE 4. Relative error on Spearman's rho 



In Table H] the relative error of the Spearman's rho is presented, when a Clayton baseline 
copula is substituted by an AMH. The error deteriorates while propagating along z. 

6. Conclusion 

In this paper a new model for bivariate dependence is proposed. The margins follow a PH 
model; the covariate acts both on the margins and on their dependence function, handled 
through the copula. Such model is in contrast with other approaches in which the role of the 
covariate is restricted to its influence on the margins, the copula deserving a separate study. 
Applications of the present model fit in "father and sons" paradigm as well as in various 
industrial contexts. Classical families of copulas (archimedean or extreme value) are shown 
to be stable under the action of the covariate, in the sens that the model acts in a transitive 
way in those classes. In parametric classes this model provides explicit transformation of the 
parameter of the copula w.r.t the covariate. This model is adequate for the case of positive 
dependence between the margins for all the values of the covariate. 



BIVARIATE COX MODEL AND COPULAS 

Appendix A. Proof of Proposition CD 
Trivially we show that lim H z (xi,x 2 ) = 0, j = 1,2; and lim H z (xi,x 2 ) = 1 
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zi-^O 



Now we will prove that H z satisfies the rectangle inequality ( Joe ( 1997f ). p. 11) which we 
recall here: for all (a 1; a 2 ), (bi, b 2 ) with a\ < bi, a 2 < b 2 , 

AH := H z (a u a 2 ) - H z ( ai , b 2 ) - H z {b u a 2 ) + H z (h, b 2 ) > 0. 



A^=(F°(a 1 ; N 



*(*) / o \ *(*)' 

H"(ai,a 2 ) ) - [H (ai,b 2 



TT(bi,a 2 )) -(Tt(bub 2 



We denote r(t) = [7T'(t, b 2 )) — ( H (t, a 2 ) ) . If ( ) is a sdf (which implies that 



H j is a 2-increasing function), then by Lemma 2.1.3 in iNelsenl (120061 ). the function r 
is nondecreasing. Therefore 



AH 



F°(6i) ) r(6 x ) - ( F°(ai ) r(a x ) 



r(6i 



r(6 



Under (JH} and if $(z) > > 0, both (iTJ and (~F°) arc sdf's. ("so the 



fact that r(&i) is negative (which holds since ( if J is a decreasing function of its second 
argument) to obtain 



> 0. 



Note that if (|H]) does not hold, # is still a sdf when $(z) > V(z) > 1. ■ 

Appendix B. Proof of Proposition [2] and [3] 
B.l. Proof of Proposition [2l Write Cjjo(u,v) = <^~ 1 (y9 (-u) + ip (v)). Using (13. 12[) some 
calculus yields (I3.15p . We prove that tp z (t) = ip (t™^M*KWT)j is also a generator for all 
values of t in [0, 1]. 

(1) Vz (l) = ¥>„(!) = 

i 

(2) <p z (i) is strictly decreasing in t, since t»i»C*W.*W) is increasing in t and (p (t) is strictly 
decreasing. 

(3) It holds ip' (t) + Vo(0 > since since C-jt is Tp 2- Let us prove that (p" z (t) > 0. We 
have, 

tp" z {t) = m{m - l)t m - 2 <p' (t m ) + {mt™- 1 ) V" (t m ) 
> m 2 t m - 2 [v' {t m ) + f>"(t m )] > 
where m = ^^nj^nyj ^ which proves that <p z is convex. 
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We now prove that C-jp in (13.15[) is TP 2 . It is readily checked that the product of two TP 2 
functions is TP 2 . Now (u,v) H- u a v b is TP 2 for all < a, b < 1 and ("U,f) i-> C iPz (u a ,v b ) is 
TP2 for all < a, b < 1, since (p z satisfies ip' z (t) + t(p"(t) > 0. We conclude that Cjj* is TP 2 
as a product of two TP 2 functions. ■ 



B.2. Proof of Proposition [3l CTg* is (I3.18|) through simple calculations. That (u,v) 1— > 
Cj[z(u,v) is TP 2 is proved as in Proposition [2] ■ 



Appendix C. Proof of Proposition H 



If K(z) < 1 then 



C H z(u,v) 



exp 



= exp 



$0) - #0) 



In w 



exp 



ln-u lnt> . , 
+ I A 



lnt; 
*(z) 



In tt I In t) 



*(*) 



In // 4- I —^-\nu + \nv } A 

$0) 



Inv 



fglnw + lnv 



exp 



(li — w s ) 



--11 i,, x ln^ + 



= exp 



exp 



lnt> 



(1 - S 



$(z) - 



$(z) 



+ 



.4 



.4 



In f 



\n.{uv)B{ 



lnv 
In 



with 



(C.27) b;(s) = (1 - + 



.4 



Hence 



= (1 - s)[l - K{z)\ + [(1 - s) J PT(z) + s] A 



(l-s)K(z) + s 



1 - K(z) - s[l - K{z)} + [(1 - s)K(z) + s}A 



(l-s)K(z) + s 
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If K(z) > 1 

C H z(u,v) 



cxp 



= exp 



(u=v s ) 



cxp 



cxp 



V(z) - $(z) 

W) 

V(z) - $>(z) 



cxp 



In u In v 



A 



In u 

*(*) 



In /• • I In u • ^Hrlu r ) 4 



*(z) 1 #(z) . 

In v 

IQlnw + lnw 



In v + 



9(z) 



i-i 

s 



= cxp 



In («■?;)-£>. 



In v 
In 



+ 



(l-s) + s 



U 1 J $(z) ^ \ 



lnt> 



with 
(C.28) 

Hence 



9(z) 



(l-s) + s 



9(z) 



A 



B 2 '(a) = 5 1 



K(z) 



1 - s 



(l-s)K(z) + s 



k ( ^- k ^ + kW) [{1 - s)K{z) + s]A 



(l-s)K(z) + s 



If = 1 

Ch*(u,v) = exp ha(uv)B^ 
with 

(C.29) 51(5) = A(s) 

We have proved that whatever K(z) 

B z (s) = 1 - mm(K(z), 1) - smin [1 - K(z)\ 



In uv 



+ min 



1 [(1 - s)K(z) + s] A 



(l-s)K(z) + s 



= 1 - W(z)K(z) - sW(z)[l - K{z)\ + W(z) [(1 - s)K(z) + s]A 

We prove that B z is a dependence function. It holds 

B z (0) = 1 -W(z)K(z) + W(z)K(z)A(0) 
= 1, since A(0) = 1 



(1 - s)tf(z) + s 



22 MOHAMED ACHIBI AND MICHEL BRONIATOWSKI 

and 

B z (l) = l-W(z)K(z) -W(z)[l- K(z)] + W(z)A(l) 
= 1 - W(z) + W(z)A(l) 
= 1, since A(l) = 1. 

We prove the upper and lower bounds for B z . 
Upper bound. Using A(s) < 1 for all s in[0, 1] 

B z (s) < 1 - W(z)K(z) - sW(z)[l - K(z)} + W(z)[(l - s)K(z) + s] = 1. 
Lower bound. Using A(s) > max(s, 1 — s) 

B z (s) > 1 - W(z)K(z) - sW{z)[\ - K(z)} + W{z) max [s, (1 - s)K(z)} . 

We prove that the RHS in the above display is larger than both s and 1 — s. 
Since max [s, (1 — s)K(z)] > s, 

RHS > l-W(z)K(z)-sW(z)[l-K(z)]+sW(z) 
= l-(l-s)min(K(z),l)>s. 

Since max [s, (1 — s)K(z)] > (1 — s)K(z), 

RHS > l-W(z)K(z)-sW(z)[l-K(z)] + (l-s)W(z)K(z) 

= 1 - sW(z) >l-8 

as sought. It remains to prove that B z (s) is a convex function. Some calculus yields 



d 2 B z 



W(z)K 2 (z) d 2 

ds 2 yS) ~ [(1 - s)K(z) + sf dt 2 



A(t) 



> 



as sought. 



Appendix D. Proof of Proposition [5] 



Write 



B z \s) = l-W{z')K{z')-sW{z')[l-K{z')]+W{z')[{l-s)K{z')+s}A 



K{z'){\ -s) + s 



In the above display it holds 



.4 



/ 



K(z')(l- s) + s 



A 



\ 



K{z) 1 
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The RHS in this latter expression can be written as a function of B z 



calculus yields 
(D.30) 

B z \s) = 1 - W{z')K{z') - sW{z')[K{z) - K{z')\ - [1 - W{z)K{z)\ 



. Some 



W(z') 
W(z) 



'l-s)^l + s\ 
K{z) 



+ 



W(z') 
W(z) 

which is (EOTD . 



* } K{z) + 



B z 



a -«)$$+«. 
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